# Identify home locations with respect FHSZ, SRA, LANDFIRE variables, and other spatial data
pacman::p_load(raster, tidyverse, fst, modelsummary, mapview, sf, qs)

source("code/globals.R")

# Load home points for at-risk sample
homepoints <- qread(file.path(WORKING, "homepoints.qs")) %>%
  select(-State) %>%
  rename(incidentid = UNIT_AND_NUMBER)

# Extract latitude and longitude so we have it downstream
lnglat <- homepoints %>% st_transform(4326) %>% st_coordinates() %>% as_tibble()
names(lnglat) <- c("lng", "lat")

homepoints <- bind_cols(homepoints, lnglat)

# Construct indicator for being in Paradise
paradise <- st_read(file.path(INPUT_PUBLIC, "ca_places")) %>% 
  filter(NAMELSAD == "Paradise town") %>% 
  transmute(in_paradise = TRUE) %>%
  st_transform(st_crs(homepoints))

homepoints <- homepoints %>%
  st_join(paradise, left = TRUE) %>%
  mutate(in_paradise = replace_na(in_paradise, FALSE))

# Assign information from WRC (see globals.R for the defn of this function)
homepoints <- assign_wildfire_hazard(homepoints)

# Assign other spatial information for all original homepoints
# A bit slow, mostly the census block part, but still runs inside of a few minutes
homepoints <- assign_spatial(homepoints)

# Set up for merging to other datasets
homepoints_df <- homepoints %>% st_drop_geometry(geometry)

# Save
write_fst(homepoints_df, file.path(WORKING, "homepoints-spatial.fst"))



# Merge original homepoints to homepoints-spatial to double check


# homepoints <- qread(file.path(WORKING, "homepoints.qs")) %>%
#   select(-State) %>%
#   rename(incidentid = UNIT_AND_NUMBER)
# test <- homepoints %>% 
#   inner_join(homepoints_spatial, by = c("ImportParcelID", "BuildingOrImprovementNumber", "incidentid"))
